#include"dividedDiff.h"

vector<vector<double>> get_diff_table(const vector<double> &inter_pts,
                                      const vector<double> &vals_of_pts)
{
    int num_of_pts = inter_pts.size();  //get the number of interpolation points

    /*定义`num_of_pts`行`num_of_pts`列的二维数组,临时记为A.实现目标是让该函数返回A*/
    vector<vector<double>> A(num_of_pts, vector<double>(num_of_pts, 0));

    /*存储差商表的第一列*/
    A[0][0] = vals_of_pts[0];
    for(int i = 1; i < num_of_pts; i++){
        /*判断传入的插值点是否相同(插值点由小到大传入),相同的插值点赋予相同的函数值*/
        if(inter_pts[i-1] != inter_pts[i])
            A[i][0] = vals_of_pts[i];
        else{
            A[i][0] = A[i-1][0];
        }
    }

    /*存储差商表的第j列, j = 1, 2,...,num_of_pts-1*/
    for(int j = 1; j < num_of_pts; j++){
        /*对于广义差商,通过指标`id`和`count`来定位传入的导数值,并存储到差商表内*/
        int id = j; //在开始计算第j阶差商时,将传入导数值定位到
                    //数组`vals_of_pts`的第j个位置
        int count = 0;  //计数:每当连续存储导数值后,count的值就+1;
                        //每当从存储导数值变为存储普通差商时,count归零.
        for(int i = 0; i < num_of_pts-j; i++){
            if(inter_pts[i+j] != inter_pts[i])  //正常差商计算
                A[i][j] = (A[i+1][j-1] - A[i][j-1]) /
                          (inter_pts[i+j] - inter_pts[i]);
            else{
                A[i][j] = vals_of_pts[id];  //存储定位到的导数值
                if(inter_pts[i+j+1] == inter_pts[i+j]){
                    id--;   //如果下一次循环仍然是计算广义差商(即计算导数值),
                            //那么由于插值点按顺序传入该次循环,则下一次
                            //要存的导数值跟这次一样,所以不必改变id的值.
                            //(这里`id--`是为了与后面的`id++`抵消)
                    count++;//下一次循环还要存储导数值,count自增.
                }
                else{
                    id += count;//如果在下一次循环计算普通差商,则将之前停留在
                                //导数值的下标id的值加回去,指向下一个正常函数值.
                    count = 0;  //导数转普通差商,归零
                }                   
            }
            id++;   //id自增,指向下一个函数(导数)值,之后再进行存储
        }
    }
    return A;
}

void print_diff_table(const vector<vector<double>> &diff_table){
    int row = diff_table.size();
    int col = diff_table[0].size();
    for(int i = 0; i < row; i++){
        for(int j = 0; j < col - i; j++){
            cout << diff_table[i][j] << "\t";
        }
        cout << endl;
    }
}

void print_diff_table_col(const vector<vector<double>> &diff_table, int col){
    int row = diff_table.size();
    for(int i = 0; i < row - col + 1; i++){
        cout << diff_table[i][col-1] << "\n";
    }
}

void print_diff_table_row(const vector<vector<double>> &diff_table, int row){
    int col = diff_table[0].size();
    for(int j = 0; j < col - row + 1; j++){
        cout << diff_table[row-1][j] << " ";
    }
    cout << endl;
}

vector<double> get_diff_table_col
        (const vector<vector<double>> &diff_table, int col)
{
    vector<double> b;
    int row = diff_table.size();
    for(int i = 0; i < row - col + 1; i++){
        b.push_back(diff_table[i][col-1]);
    }
    return b;
}